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NOMENCLATURE 

ollov/ing  nomenclature  is  used  in  the  paper: 

Velocity  of  sound  in  fluid 

Joukowsky  water  hammer  velocity  (See  eq,  l) 

Bulk  modulus  of  fluid 

Radius  of  tube 

Thickness  of  tube  wall 

Young’s  modulus  of  tube  wall  material 

Pressure  in  fluid 

Axial  velocity  of  fluid 

Distance  axially  along  the  tube 

Time 

Axial  displacement  of  tube 

Radial  deflection  of  tube  wall  (+  outward) 

Initial  density  of  fluid 

Radial  distance 

Velocity  potential 

Mass  of  tube  per  unit  surface  area 

Poisson’s  ratio  for  tube  wall 

Velocity  of  sound  in  tube  wall 

Force  per  unit  area  on  tube  wall 
Transform  variables 
Integers 
See  eq.  (23) 

See  eq.  (2li) 

Bessel  functions 

See  eq.  (30) 

See  eq,  (36) 

See  eq.  (37) 

See  eq,  (38) 


C',C 

Arbitrary  constant  velocities 

I,  H 

Integrals,  See  eq. 

°1»°2 

Phase  velocities 

L 

Length  of  wave  front 
c  l 

*• 

See  eq,  (U8) 

R  = 

0 

Variable  of  integration 

c 

££ 

m 

p 

See  eq.  ($0 

A  = 

1.  Introduction 

The  usual  theory  of  water  hammer  which  is  used  to  compute  the  magnitude  and 
velocity  of  pressure  waves  in  a  pipe  filled  with  an  elastic  fluid  was  given  by 
Joukowsky  [8]  in  1898.  This  theory  predicts  that  pressure  waves  travel  without 

change  of  shape  at  a  velocity  c  s 
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C 


e  " 


C 


_l  2K a 

"TfT 


(l) 


In  the  derivation  of  this  formula  it  is  assumed  that  the  pressure  is  uniform 
across  any  section  of  the  pipe  and  that  the  deflection  of  the  pipe  is  equal  to  the 
static  deflection  due  to  the  instantaneous  pressure  in  the  fluid.  The  first 
assumption  neglects,  in  effect,  the  inertial  forces  associated  with  radial  motion 
of  the  fluid.  The  second  neglects  the  mass  of  the  pipe  wall  and  the  longitudinal 
and  bending  stresses  in  the  pipe  wall. 


For  the  case  of  an  infinite  train  of  sinusoidal  pressure  waves  in  a  tube  filled 
with  fluid,  much  work  has  been  done  without  resorting  to  the  above  assumptions.  In 
1878  Korteweg  [9]  showed  that  the  phase  velocity  of  sinusoidal  waves  varies  with 
the  wave-length  when  the  radial  inertia  of  the  fluid  and  the  mass  of  the  pipe  wall 
are  taken  into  account,  Korteweg  also  gave  the  expression  for  c  ,  eq.  (1)  as  an 
approximation  of  the  phase  velocity  for  long  wave-lengths.  The  work  of  Korteweg 
was  extended  in  1898  by  Lamb  [10]  who  considered  the  pipe  as  an  elastic  membrane 
thus  including  the  effect  of  longitudinal  stresses  in  the  pipe  wall.  Lamb  showed 
there  are  two  finite  phase  velocities  as  the  wave  length  approaches  infinity. 

One  of  these  is  very  close  to  Korteweg' s  result  c  and  the  other  is  very  nearly  the 
velocity  of  sound  in  the  material  of  the  tube  wall.  Recent  papers  by  Jacobi  [6] 
and  Thomson  [lh]  include  the  flexural  stiffness  of  the  pipe  wall.  The  inclusion 
of  the  bending  stresses  affects  the  phase  velocities  appreciably  only  for  very  short 
wave-lengths,  i 

in  .in cketP  refer  to  the  Bibliography  at  the  end 


In  the  present  paper  water  hammer  waves  are  considered  without  some  of  the 
simplifications  made  in  the  Joukowsky  theory.  The  equations  of  motion  for  a 
thin  cylindrical  tube  are  used,  cf.  Fliigge  [11],  ^ 

with  one  additional  term  for  the  effect  of  rotatory  inertia.  The  fluid  is  assumed 
to  be  elastic  and  inviscid  but  no  assumptions  are  made  regarding  its  motion  except 
that  it  is  small. 

Only  the  case  where  the  velocity  of  sound  in  the  material  of  the  pipe  wall  is 
greater  than  in  the  fluid  will  be  considered.  Much  of  the  analysis  would  be  valid 
in  the  opposite  case,  but  the  relative  speeds  of  the  various  waves  would  differ* 

2.  A  Formal  Solution  by  Integral  Transforms 

obtained  in  the  form  of  a.  double  integral. 
In  this  section  a  formal  solution  of  a  typical  water  haoner  problem  is^ 

The  solution  is  reduced  to  terms  of  single  real  integrals  whose  evaluation  is 
discussed  in  the  next  section. 

The  problem  considered  is  to  find  the  motion  of  a  thin  cylindrical  tube  filled 
with  fluid  starting  from  given  initial  conditions.  The  conditions,  shown  in  Fig.  1, 
are  that  for  t  ■  0: 

p  ■  0,  v  -  0  for  z  >  0  'j 

P  »  P0,  v  -  vo  for  z  <  0  >  (2) 

u  ■  0,  w  ■  0,  for  all  z  J 

The  velocity  vq  is  assumed  to  be  related  to  the  pressure  pQ  by: 


This  is  the  velocity  that  would  exist  behind  a  plane  shock  wave  of  pressure  po  in 
an  infinite  body  of  fluid.  Thus  the  initial  conditions  for  the  fluid  correspond  to 
a  step  shock  wave  moving  in  the  positive  direction.  An  external  pressure  equal  to 
PQ  is  assumed  acting  on  the  outside  of  the  tube  for  z  <  0  to  prevent  any  motion  of 
the  tube  at  z  ■  -  oo  .  This  external  pressure  is  assumed  to  remain  on  the  tube  for 
z  <  0  for  all  t. 


These  initial  conditions  are  intended  to  simulate  the  situation  arising  when 
a  pipe  leading  from  a  reservoir  is  filled  with  fluid  at  rest  below  reservoir 
pressure  and  a  valve  at  the  reservoir  end  is  suddenly  opened.  The  part  of  the  tube 


left  of  the  origin  in  Fig,  1  is  a  convenient  device  for  supplying  fluid  and  energy 
in  place  of  a  reservoir.  It  is  expected  that  the  particular  means  of  generating  a 
water  hammer  wave  will  not  essentially  affect  its  characteristics  when  the  wave  has 
progressed  far  enough  down  the  pipe. 


The  motion  of  the  fluid  will  be  described  by  a  velocity  potential  which  must 
satisfy  the  equation: 

CM 


v24>  -  <b  +  =  -L-X 

T  Yrr  r  Yr  C2-  Y 


where  the  subscripts  r  and  z  denote  derivatives. 

The  motion  of  the  tube  due  to  a  radial  force  per  unit  area,  q  (z,  t),  is 
governed  by  two  equations: 

a  '  ■  ‘  ‘ 
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o 

detailed 

A  ^derivation  of  these  equations  is  given  in  [12],  They  agree  with  the  equations 
given  by  Fliigge  [11]  except  that  the  rotatory  inertia  term  w__  has  been  added  above. 

The  fact  that  the  fluid  and  tube  are  expected  to  remain  in  contact  is  expressed 
by: 


w  =  — 


dr 


Jr*  a 


(7) 


The  pressure  exerted  on  the  tube  wall  by  the  fluid  is: 


r-a 


(8) 


This  is  the  value  of  q  to  be  used  for  z  >  0  in  eq.  (*>).  For  z  <  0,  the  external 

pressure  -  p  must  be  added, 
o 


It  is  convenient  to  consider  the  solution  of  the  problem  as  the  sum  of  two  parts 


The  first  part  is  the  shock  wave  that  would  be  produced  by  the  initial 
conditions  if  the  tube  were  rigid.  The  tube  displacements  Wp  Up  for  this  part 
are  zero.  The  velocity  and  pressure  Vp  pp  for  the  fluid  are: 


vt=v0 

PrPo 


for  2  *ct 


for  2  >  ct 


(9) 


The  second  part  is  the  response  of  the  tube  and  fluid  to  a  progressive  load 
pressure  p  equal  to  the  pressure  produced  by  the  shock  wave  of  the  first  part.  This 
load  pressure  acts  on  the  tube  wall  radially  outward.  Its  magnitude  is  p  for 
0<  z  <  ct  and  zero  elsewhere. 


Inasmuch  as  Wp  Up  Vp  etc.  are  known,  the  problem  reduces  to  finding  the 

the  response  w2,  u2,  Vp  etc.  of  the  tube  and  fluid  to  the  load  pressure  p8.  For 
this  second  part  of  the  solution  the  system  is  Initially  at  rest  and  the  governing 
equations  are  (1*),  ($),  (6),  (7)  with: 

+  Ps  (10) 

r=a 

The  solution  is  begun  by  transforming  the  governing  equations  and  eliminating 
the  velocity  potential  <f>2«  The  remaining  transformed  equations  are  solved  simultane¬ 
ously  and  the  inversions  of  the  transformations  then  yield  the  solution. 


Po 


3  4> 

3t 


With  respect  to  z  a  Fourier  transform,  cf.  [15]  p.  i*2,  is  used.  Defining  the 
transform  of  a  function  g  (z,  t)  with  respect  to  z  by: 


(ii) 


then  the  inversion  is: 


■  CD 


OO 


(12) 


J 

-OO 


With  respect  to  time  a  Laplace  transform  is  used,  but  the  transform  variable 
adopted  is  -  L  -  -J~~  1  times  the  usual  transform  variable.  This  results  in  a  certain 
symmetry  with  the  first  transform  and  convenience  in  later  analysis.  Defining  the 


(13) 


transform  of  a  function  g  ( ^  ,  t)  with  respect  to  t  by: 

OO 


dt 


then  the  inversion  is: 


OO-lot 


Cli»> 


-co-  cok 


The  double  transform  of  derivatives  of  g  (z,  t)  may  be  expressed  in  terms  of  the 
double  transform  of  g  (z,  t)  by  the  usual  partial  integration  procedure*  [133  P*  27. 
Using  the  fact  that  initial  velocities  and  displacements  are  zero*  the  transforms  of 
the  derivatives  involved  are: 


y 

a*"  dts 


=  (^)n(iiOs 


(15) 


Using  q2  as  given  by  eq,  (10)*  the  double  transformation  of  eqs,  (U),  (5)*  (6) 


and  (7)  yields: 
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(16) 


(17) 


(18) 


(19) 


The  double  transforms  are  functions  of  \  and  XL  only  except  for 
$  2  which  is  4>2  (  f  ,  XL,  r).  From  eq,  (16)  it  is  found  that  the  r  dependence  of  $  2 
involves  a  Bessel  function  of  the  first  kind.  The  function  pg  is  known.  Solving 
eqs.  (16),  (17),  (18),  (19  for  *2,  $2  and  inverting  yields: 


oo 


00-i.ok 


The  evaluation  of  the  radial  deflection  will  be  considered  first*  It  is 
typical  of  the  quantities  of  interest#  The  inner  or  XL  integral  of  eq,  (20)  nay  be 
evaluated  by  Cauchy's  residue  theorem  applied  in  the  complex  Xu  plane#  During  this 
process  \  is  regarded  as  a  constant  real  parameter  since  the  outer  Integral  involves 
only  real  values  of  ^  . 


The  factor  XL  ^  in  the  denominator  of  the  Integrand  in  eq,  (20)  produces  a  simple 

pole  at  -Cl  -  0  because  (  XL2  I)  is  finite  at  XL  -  0.  There  is  no  pole  at  (&♦  f  )  -  0 

N  Ci 

because  (~  ♦  \  )  A  is  finite  at  this  point.  There  are  poles  at  the  points  for  which 
the  remaining  bracket  is  zero,  i,e,,  where: 


O 


(25) 


These  roots 


For  any  real  \  this  equation  haa  an  infinite  number  of  roots  n« 
have  an  interesting  physical  significance.  They  are  the  circular  frequencies  of  the 

inodes  of  free  vibration  of  the  tube  and  fluid  system  having  a  wave-length  A  ■  4"  • 

£l-  ' 

Alternatively  the  ratios  are  the  phase  velocities  cn  of  progressive  sinusoidal 

These  facts  may  be  developed  by  considering  the  governing 


2tt 


waves  of  wave-length  y  . 

equations  (U),  (5),  (6)  and  (7)  with  q  in  eq.  (5)  equal  to  the  fluid  pressure  only. 
A  solution  is  assumed  of  the  form: 


W 
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<t>  = 


,  e  L  ( f  2-  +nt) 
,  c  i{<-i  +nt) 

if  f(r)  e 


w 
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(26) 


where  the  prime  quantities  are  constants  and  f  (r)  is  an  unknown  function  of  r. 
Substituting  eqs.  (26)  in  eqs.  (U),  (55,  (6)  and  (7)  the  constants  and  f(r)  may  be 
eliminated,  leading  to  exactly  the  transcendental  eq.  (25)  relating  ^  and  fL  . 

For  real  values  of  f  all  the  roots  fL  Q  of  eq.  (25)  are  also  real.  If  any  of 
the  XL  n  were  complex  or  imaginary  the  corresponding  free  vibration  would  grew  or  decay 
indefinitely  in  time.  Since  the  system  has  no  means  of  dissipating  or  radiating  energy, 
this  is  not  possible  physically. 

Fig.  2  shows  the  XL  vs  j  values  for  a  typical  case  in  which  the  velocity  of  sound 
in  the  tube  wall  is  greater  than  in  the  fluid,  i.e.,  c  >  It  is  found  that  for  any 
value  of  |  there  is  one  root  XI,,  for  which  (*=j  -  y )  is  negative.  This  root  corre¬ 
sponds  to  a  phase  velocity  less  than  c,  the  velocity  of  sound  in  the  fluid.  There  are 

q2  2 

also  an  infinite  number  of  roots  for  which  (-=^  -  ^  )  is  positive.  These  correspond 

c 

to  phase  velocities  greater  than  c.  Fig.  3  shows  the  same  data  as  in  Fig.  2  plotted  as 
phase  velocity  vs  wave-length. 

The  fact  that  two  sets  of  XL  curves  pass  through  the  origin  in  Fig,  2  shows  that 
there  are  two  modes  which  have  finite  phase  velocities  as  the  wave  length  approaches 
infinity  and  the  frequency  approaches  zero.  This  was  also  shown  by  Lamb  [10],  p,  7. 

It  is  mentioned  here  because  in  a  recent  paper  [Uj],  p.  933,  it  is  erroneously  stated 
that  there  is  a  cutoff  frequency  below  which  only  one  mode  may  be  propagated. 


The  above  discussion  shows  that  the  poles  of  the  integrand  of  eq.  (20)  are  at 

XL  •  0  and  XL  -  xl.  .  all  on  the  r. al  xl  axis.  The  path  of  integration  is  below 

n 

this  axis  and  it  will  be  closed  in  the  usual  manner  by  adding  a  semicircle  of 
infinite  radius  in  the  upper  half  of  the  XL  plane.  The  integral  over  the  semicircle 

approaches  zero  as  the  radius  approaches  infinity.  The  contribution  of  the  pole  at 

St  •  0  to  w2  is: 


Oo 


The  contributions  of  the  Xl^  poles  to  w^  may  be  expressed  as: 

00 


The  sum  indicated  is  over  the  XL  roots  for  each  value  of  f  ,  The  values  of  the 

n  > 

as  functions  of  \  are  given  by  Fig,  2,  Eq.  (28)  may  be  written  in  the  form 


The  integral  in  eq.  (29)  covers  all  four  quadrants  of  the  XL  curves  in  Fig,  2, 
but  it  may  be  expressed  as  an  integral  over  the  first  quadrant  values  only  because 
D  is  even  in  Xl  and  even  in  ^  and  the  xl  curves  are  symmetric  about  both  axes. 

Thus: 

00  _ 

stn(yt-XLt) _ s vniU  +cl±)  .  , 

(t-flD  "  (fH)D  J  5  (31) 

0  Xl  n 

where  the  summation  is  over  the  XL  curves  in  only  tne  first  quadrant  of  Fig,  2,  This 
concludes  the  formal  solution  for  w.  The  evaluation  of  the  integrals  remaining  in 
eqs,  (27)and  (31)  will  be  considered  in  the  next  section. 

Instead  of  evaluating  and  §  certain  quantities  of  interest  which  are 


10 


derivatives  of  u2  and  $  2  will  be  computed  directly.  These  are: 

Ha  ^ 

Pressure  In  the  fluid:  Pa=V»Tt” 

^$>z 

Axial  velocity  of  the  fluid: 


dt 


(32) 


Longitudinal  strain  in  the  tube: 


u«sir 


Expressions  for  p2#  v2  and  uz2  are  derived  from  the  double  integrals  for  u2  and  <|>  2, 

eqs.  (21)  and  (22),  through  application  of  eq,  (15).  The  inner  integral  in  ft  may¬ 
be  evaluated  in  each  case  by  using  Cauchy’s  residue  theorem  as  in  the  case  of  w2# 

The  results  are: 
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-(t-f)B  (-f +f)B 

n  . 


d.\  (35) 


(36) 


(37) 


(38) 


and  f(|)  stands  for  the  denominator  of  the  integrand  in  eq,  (27) 
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The  solution  of  the  problem  originally  posed  at  the  beginning  of  this  section 
is  given  formally  by: 


w  “  w2a  *  w2b 


P  "  Pi  ♦  Pi 


v  “  V1  ♦  y2 


u 


u 


z2 


Eqs.  (27)  and  (31) 
Eqs.  (9)  and  (33) 
Eqs.  (9)  and  (3U) 
Eq.  (3$) 


(39) 


3.  Evaluation  of  the  Solution  for  Large  Values  of  |z[  and  t. 

In  the  formal  solution  indicated  by  eqs.  (39)  the  integrations  with  respect  to  | 
would  be  difficult  to  carry  out  in  general  because  the  XL  values  are  known  as 
functions  of  ^  only  through  the  transcendental  equation  (25).  However,  some  informa¬ 
tion  may  be  obtained  by  considering  approximations  of  the  integrals  for  large  values 
of  | z  |  and  t. 


Consider  the  part  of  w  given  by  eq.  (31)  and,  in  particular,  values  of  a  and  t 
such  that: 


*  «  c'  t 

'die re  c*  is  an  arbitrarily  selected  velocity.  Substituting  in  eq,  (31)* 

00 

W  =  — — 

2b  IT  m 


(1*0) 


(y 

sm  C 

-XL^t 

Sin  (fc1 

L_ 

K 

(T-f 

)D 

/  fL 
(—  + 

. 

(la) 


It  is  found  by  inspection  of  D,  that  the  functions 


and 


c  i  '  ~  '  c  f  ^  D 

are  bounded  for  all  of  the  XL  curves  for  all  \  except  for  the  lowest  two  curves,  XL, 
and  XL  2»  the  point  ^  *  0,  The  distinction  between  the  lowest  two  and  the 
remaining  XL  n  curves  arises  because  only  the  XL  ^  and  XL  g  curves  pass  through  the 
origin.  A  physical  interpretation  of  this  fact  is  that  only  these  two  lowest  modes 
have  finite  phase  velocities  as  the  wave-length  increases  indefinitely. 


For  the  higher  XL  n  branches,  n  >  2,  the  integral  (la)  may  be  approximated  by 
the  method  of  stationary  phase  [7]  p.  505,  for  large  values  of  t.  This  method  shows 


12 


these  integrals  to  be  of  the  order  of  — —  • 

sft 


The  integrals  over  the  XL  ^  and  _Tl  ^  curves  will  be  broken  into  two  parts;  one 
part  from  |  =  0  to  [  =  6  and  a  second  part  from  €  to  oo  ,  where  6  is  a  small 
positive  number.  By  the  stationary  phase  theorem  the  portions  from  6  to  oo  are  again 
0  (  — -  ).  Collecting  these  results,  eq.  (1*1)  becomes: 


/ T 


Po  f 


irm 


I 


XL 


sinftc'-iOt  sin 


(t-O*  (t+0D 


dl  +0(/r) 


0*2) 


Since  €  may  be  chosen  indefinitely  small,  the  functions  in  the  integrand  may 
be  approximated  by  an  appropriate  number  of  terms  of  their  expansions  about  ^  "  0. 
Consider  the  first  term  of  the  integral  with  XL  »  XL  1  which  may  be  written: 


w  =  A 
I  irm 


Sin  CfC-lpt 

!  Kf-p 


d! 


(1:3) 


The  denominator  of  the  integrand  approaches  a  constant  times  ^  as  j  approaches 
zero.  The  constant  is: 


where: 


f- 


/_£i_  \ 

)D> 

l  c  1 

C,  =  lim  (  ,  ) 
*-*-0  N  1  ' 

+e0ct 

*c02 

(-# 

0*1*) 


\ 


4-1  -  V 


‘+Jl 


12  <r 


\  0*5) 


J 


The  constant  c^  is  the  phase  velocity  of  the  slowest  waves  for  wave-lengths 
aching  infinity*  It  is  very  < 
for  c^  is  derived  in  the  appendix* 


approaching  infinity*  It  is  very  close  to  the  Joukowsky  velocity  c  .  An  expression 

© 


In  the  numerator  of  eq.  (1*3)  XL^  is  replaced  by  the  first  two  tenas  of  its 


expansion  about  ^  "0: 


(46) 


This  expansion  is  developed  in  the  appendix.  It  is  necessary  to  use  two  terns 
of  the  expansion  because  the  first  tern  alone  would  result  in  a  vanishing  numerator 
for  the  case  c*  ■  c^. 

Using  eqs.  (UU)  and  (1*6)  in  (1*3)  yields: 

00 

W  = _ B> _ 

0 

where  z‘  is  defined  by: 

^(C-C^t  (U8) 


(1*7) 


The  new  coordinate  2*  is  the  location  of  any  point  relative  to  an  origin  moving  at 

velocity  c1 .  In  eq.  (1*7)  the  upper  limit  has  been  changed  from  6  to  oo  .  The 

integral  thus  added  is  0(— by  the  stationary  phase  theorem,  and  is  negligible 

St 

in  the  present  approximation. 


For  z1  *  0,  the  integral  in  eq.  (1*7)  has  a  constant  value  of  £  provided  d^  is 
positive  which  is  the  case  here.  For  points  where  | z 1  {  is  large,  the  term 
containing  ^  **  may  be  neglected.  Then  the  integral  in  eq,  (1*7)  becomes  -  ^  for  -  z*. 
For  points  near  2*  »  0  it  is  convenient  to  change  the  variable  of  integration  from  | 
to  rj  •  i  s’.  Then: 

•  CO 

Sin  (17+  0t73)  ^ 

- — dr}  (i !*9 ) 


wi  = 


AT  m 


(■r  -OD,  J 


where: 


d,  t 


(*) 


Then: 


From  eq.  (1*9)  it  may  be  seen  that  for  w^  constant  (J  must  be  constant,  say  ^  Q. 


(51) 


This  shows  that  points  for  which  Wj  is  constant  continually  move  away  from  the  point 
2»  •  0.  To  summarize,  Wj  is  a  wave  having  constant  deflection  at  the  point  z*  »  0 
i.e.,  moving  with  velocity  c^.  At  points  far  ahead  and  far  behind  this  wave  front 
point  the  deflection  is  also  constant,  but  the  wave  front  continually  spreads  out 


around  the  point  as'  *  0.  Fig,  4  shows  the  general  shape  of  the  wave. 


If  the  same  reasoning  that  has  been  used  to  evaluate  Wj  is  applied  to  the 
remaining  terms  of  eq,  (1*2)  similar  results  are  obtained.  Thus  for  large  t: 


W2b=  II  "  Hl  +  I2  " 


where: 


I, 

H, 


irm 


L 


n 


Sm(i  +  CMt  t)f 

\ 


(53) 


The  constants  c 2,  Dg,  d2  v’ich  appear  in  H2  and  1^  are  defined  in  the  same 
way  as  c^,  D^,  in  eqs,  (1*5)  and  (1*6)  by  changing  all  subscripts  from  1  to  2  in 
these  equations.  The  velocity  c2  is  very  nearly  the  velocity  of  sound  in  the  tube 
wall,  cQ. 


The  term  1^  is  the  deflection  Vj  discussed  above.  The  term  is  a  wave  of 
the  een.  velocity  end  .hep.  a.  .j  hut  of  —Her  edited,  end  traveling  in  the 
negative  direction. 


The  H2  and  I2  terms  are  also  waves  of  the  same  shape  as  w^,  but  traveling  at 
the  velocity  c2«  These  precursor  waves  are  found  to  be  very  small  compared  to  the 
main  water  hammer  waves  1^,  which  travel  at  velocity  c^. 


The  complete  deflection  w  is  the  sum  of  w^,  eq.  (52)  and  w2a,  eq.  (27),  The 
deflection  w2ft  is  independent  of  time  and  its  value  for  large  values  of  |zj  may  be 
approximated  by  the  same  kind  of  reasoning  as  used  for  w^  for  large  values  of  t. 
The  result  is: 


W 


_  + 


2<X 


2  E  K 


I  + 


II  <L 


(\-vl) 


for  i  * 


(51*) 


The  term 


12a2 (1  -  lA) 


walled  tubes  are  considered 


will  be  neglected  compared  to  unity  since  only  thin- 
•  Thus: 


for  ±  2 


-  If- 


m! 

2Eh 


m 


The  procedures  used  to  evaluate  v  may  also  be  used  to  obtain  approximations  of 
p,  v  and  u  for  large  |zj  and  t  starting  with  the  general  expressions,  eqs.  (39)» 


The  results  are: 


W  =  t 


p= 


V  = 


Po_ 

2 


*e.c 


p  * 

—Z—  +  1 
2E  h  ' 

2PoC^ 


r  rl 


a 


Lc 

2C, 

a 

O 

-%  -1 

.  c* 

poira 

u  =  ;  — - + 

*  +  2Eh 


_C^ 


(1,+H.V 


Hz 

■for 

2P»cz 

a 

r  r*  i 
■—  -l 

c 

2C, 


44 -»i 

Lc*  J 


(V»,) 


>  (#) 


7^ 


Lc*  J 


(i*-ho  wn 


These  equations  lead  to  four  wave  fronts  in  each  case  as  discussed  for  v. 
Numerical  restats  for  a  typical  case  of  water  in  a  thin  steel  tube  are  shown  in 
Fig*  5*  In  the  present  approximation  for  large  t,  the  pressure  and  fluid  velocity 
are  again  independent  of  r,  i,e,,  they  are  constant  across  any  section* 

In  order  to  estimate  the  length  over  which  a  given  wave  front  extends,  the 
slope  w2  will  be  computed  at  the  points  which  move  with  constant  velocity  c^  or  Cg* 
The  difference  in  w  before  and  after  a  wave  front  divided  by  the  slope  at  the  wave 
front  point  will  be  indicative  of  the  length  of  a  wave  front* 


The  double  integral  for  wz  is  derived  from  w,  eq.  (20).  The  first  integral  may 
be  evaluated  and  the  result  folded  by  the  same  proced_-es  used  for  w*  For  a  point 
moving  so  that  z  »  c^  t: 


00 


o 


Po 

(r 

^cos(c,^  -Sl)t 

irm  : 

0 

L 

(tM)b 

^C0S(C^-Kfl)t 

(M)D 
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The  function  ~J^y  is  bounded  for  all  \  and  the  first  integral  in  eq.  (57) 
is  0  (^)  by  the  Itiemann-Lebesgue  lemmas  [16]  p.  172. 


In  the  second  integral  the  integrand  is  also  bounded.  This  Integral  may  be 
approximated  by  the  method  of  stationary  phase.  The  dominant  part  of  the  integral 
will  come  from  the  vicinity-  of  stationary  phase  points  which  are  defined  by* 


(|c,±n)  =  o  L.e. 


an 


=  *c. 


(58) 


In  general  such  points  give  a  contribution  which  is  0  ( — ),  but  if* 

✓  t 


•(^c,  ±n)  =  o  i.e. 


_  o 

J  c'i-  u 


(55) 


then  the  contribution  to  the  Integral  from  the  vicinity  of  this  point  is  0  ( j-—  ). 
[53  p.  1*0.  Eqs.  (58)  and  (59)  are  both  satisfied  on  the  XL  ^  curve  at  \  *  0,  * 
{fence  this  contribution  will  predominate  over  all  others  for  large  values  of  t. 


W*= 


irm('c‘~l)I)i 

This  Integral  (61)  is  known,  [5]  p.  1*1.  Hence* 

P.  r  (|)  sun  $ 


cos  i-jMij 


wt  = 


3,Tm  ('T’~  *) 


(») 


To  evaluate  the  dominant  term,  the  expression  for  eq.  (1*6) ,  is  substituted  In 
eq.  (57).  Hence  for  large  t*  „ 

I  (M*. 

6 

Changing  the  variable  of  integration  to  ^  d^  t  j  : 

oo 

Po 


(61) 


(62) 


Now  the  change  of  w  from  a  location  far  behind  the  point  *  *  c.  t  to  a  location 

of  th€  point 

far  ahead  z  •  c,  t  is  by  eq.  (1*9): 
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&W  — 


(63) 


Defining  the  "nominal  length"  1^  of  the  wave  front  to  be  aw  divided  by  w^: 


(6U) 


Computation  of  the  "nominal  length",  Lg,  of  the  wave  fronts  moving  at  velocity  Og 
results  in  the  same  formula  with  replaced  by  dg. 


Some  numerical  values  of  and  Lg  are  given  in  Table  1.  It  may  be  noted  that 
although  the  length  of  the  wave  front  is  long  compared  to  the  diameter  of  the  pipe 
it  is  short  compared  to  the  distance  from  the  origin  that  the  wave  has  traveled* 


The  above  computations  of  the  lengths  of  wave  fronts  hold  also  for  p,  v  and 
u£  because  they  are  described  by  the  same  integrals,  1^,  H^,  lg,  Hg. 

U*  An  Extension  of  Joukowsky's  Method 


Some  of  the  results  derived  by  use  of  integral  transforms  may  also  be  derived 
by  a  much  simpler  procedure  which  is  an  extension  of  Joukowsky's  method.  The 
extension  consists  of  taking  into  account  the  longitudinal  stresses  and  longitudinal 
inertia  of  the  pipe  wall*  This  makes  possible  the  precursor  type  wavr  which  is 
essentially  a  tension  wave  in  the  pipe  wall* 

The  following  derivation  follows  the  general  plan  of  Joukowsky’s  paper  [8], 

The  pressure  and  axial  velocity  of  the  fluid  are  assumed  to  be  uniform  over  any 
cross-section  of  the  pipe*  The  fluid  is  assumsd  compressible,  its  bulk  modulus 
being  K.  The  equation  of  continuity  may  then  be  written  in  the  form* 

_cW  _  J_  _bp  _2_  dw 
K  bt  *  a  bl 


(65) 
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The  equations  of  motion  of  the  fluid  reduce  to: 

,iv=__SP  (66) 

dt  bi. 


So  far,  these  are  the  equations  usually  used.  The  new  element  is  the  use  of  the 
equations  of  motion  of  the  tube  considered  as  an  elastic  membrame  instead  of  a  ring 
of  no  mass.  These  equations  are  derived  from  eqs,  (5)  and  (6)  by  dropping  all  terms 
arising  from  bending  stiffness  and  rotatory  inertia  and  substituting  p  for  q: 


P 

mc£ 


=  0 


(67) 

(68) 


In  eq,  (67)  the  term  w  has  been  omitted  because  the  effect  of  radial  inertia  of  the 
tube  and  of  the  fluid  is  neglected  as  in  the  Joukowsky  theory* 


No  attempt  will  be  made  to  solve  the  equations  (65),  (66),  (67)  and  (68)  in 
general.  It  is  sufficient  for  the  purpose  at  hand  to  show  that  the  system  permits 
solutions  which  are  waves  of  arbitrary  shape  moving  at  either  of  the  velocities  c^  or 
Cg  without  dispersion*  A  solution  is  asstaned  of  the  form: 


W  =  W1  f  (i  +  C  X)  N 
P  =  p'  f  (l  +  C  t) 

V  =  V*  f  (l  +  Ct)  f 

ui=  + 

Ut«  j 


(69) 


where  the  primed  quantities  are  constants  and  c  is  an  unknown  constant  velocity.  The 
constants  u£,  u|  and  c  are  related  by  the  requirement: 


(70) 


Substituting  the  assumed  form  (69)  into  eqs.  (65),  (66),  (67),  (68)  and  (70),  there 
result  five  linear  homogeneous  equations  involving  the  five  primed  constants.  In  order 
that  these  equations  have  a  solution  different  from  zero  it  is  necessary  for  the 
determinant  of  the  coefficients  to  be  zero*  This  condition  yields  an  equation 


_  .  _  +  ♦ 

involving  c  which  can  be  satisfied  only  if  c  takes  on  the  values  -  c^  or  -  Cg,  Thus 

the  velocities  at  which  waves  may  travel  without  dispersion  in  the  present  considera¬ 
tion  are  the  same  as  the  principal  velocities  derived  in  the  integral  transform 
method* 


A  further  point  of  interest  is  that  the  ratio  of  the  change  in  w  to  the  change 
in  p  across  a  wave  front  is  the  same  in  the  extended  Joukowsky  theory  and  the 
transform  solution.  The  same  is  true  of  any  two  of  the  quantities  w,  p,  v,  u^,  etc. 
This  is  shown  by  solving  for  the  ratios  of  interest  in  each  theory  with  c  m  c^  or  Cg 
as  the  case  may  be* 

These  results  are  of  interest  because  they  indicate  the  significance  of  the 

assumptions  made  in  the  Joukowsky  theory.  It  shows  that  c.  differs  from  c  because 

longitudinal  of  the  pipe.  e 

the  Joukowsky  theory  neglects  longitudinal  stresses  and  inertia.  It  also  shows 
that  the  dispersion  predicted  by  the  transform  method  is  due  to  the  effects  of  the 
bending  stiffness  and  of  the  radial  inertia  of  the  pipe  and  fluid*  The  effect  of 
radial  inertia  is  much  more  important.  After  the  wave  disperses  somewhat,  bending 
effects  are  negligible. 


5.  Other  Initial  Conditions 

The  integral  transform  method  may  be  applied  to  certain  other  initial  conditions 
beside  those  assumed  in  Sections  2  and  3*  The  problem  and  results  shown  in  Fig,  6 
illustrate  a  second  example.  The  initial  conditions  are: 

v  ■  v  for  z  <0  vs-v  for  z  >0 

o  o 

p»u*w=ut*w^,»0  for  all  z. 

Under  these  conditions  each  half  of  the  pipe  simulates  the  water  hammer  problem  in 
which  an  established  flow  is  suddenly  cut  off  by  an  instantaneous  valve  closure* 

Since  by  symmetry  no  fluid  crosses  the  plane  z  ■  0,  each  half  of  the  pipe  acts  as 
a  closed  valve  as  far  as  the  other  half  is  concerned. 

The  magnitudes  and  velocities  of  the  waves  may  also  be  established  on  the  basis 
of  the  simplified  procedure  of  Section  In  general,  this  will  be  the  simpler 


method*  It  does  not  predict  dispersion,  but  it  may  be  assumed  that  the  dispersion 
characteristics  are  the  same  in  any  case* 

6.  Conclusions 

The  theory  developed  herein  shows  the  following  results  which  are  not  contained 
in  the  usual  theory: 

1*  If  the  initial  conditions  specify  a  sharp  wave  front  of  pressure,  this  front 
will  gradually  disperse  over  a  finite  length  as  it  proceeds.  This 
dispersion  continues  indefinitely  so  that  no  steady  wave  shape  is  ever 
reached* 

2*  Two  waves  will  be  generated  in  general;  one  which  is  similar  to  the 

Joukowsky  result  and  one  which  is  essentially  a  tension  wave  in  the  pipe 
wall.  The  pressure  change  due  to  this  precursor  is  very  small* 

a 

3*  AfterAsufficiently  long  time  has  elapsed  since  the  generation  of  a  wave 
there  will  be  one  point  in  each  wave  front  which  has  a  constant  deflection 
and  fluid  pressure  and  moves  at  a  constant  velocity.  This  velocity  is 
very  close  to  the  Joukowsky  velocity  for  the  main  pressure  wave*  For  the 
precursor  this  velocity  is  slightly  less  than  the  velocity  of  sound  in 
the  pipe  wall. 

U.  The  bending  stresses  in  the  pipe  wall  have  very  little  effect  on  the 
characteristic  velocities  and  the  rates  of  dispersion  after  a  sufficient 
time  has  elapsed.  Such  stresses  can  only  be  significant  in  the  early 
stages  of  the  motion  in  regions  where  the  wave  front  is  sharp  so  that  the 
variation  of  pressure  and  deflection  is  very  rapid  along  the  length  of 
the  pipe. 

5*  The  relations  between  fluid  velocity,  pressure  and  hoop  stress  given  by  the 
Joukowsky  theory  are  substantially  correct,  but  there  are  also  longitudinal 
stresses  in  the  pipe  not  predicted  by  the  usual  theory. 


1 

X 


U] 

6,  Available  experimental  data  is  in  very  close  agreement  with  the  predictions 
of  the  Joukowsky  theory  as  to  the  pressure  rise  and  velocity  of  the  main 
pressure  wave.  Inasmuch  as  the  analysis  presented  in  this  paper  gives 
essentially  the  same  results  as  the  Joukowsky  theory  for  these  items,  it  is 
confirmed  experimentally.  However,  there  is  no  experimental  data  available 
on  the  precursor  wave  or  the  rate  of  dispersion  predicted  herein. 
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APPENDIX 


Expansion  of  il  ^  and  XL  2  About  ^ •  0 

The  values  _ru  XL  2  are  defined  as  the  lowest  two  roots  of  eq.  (25).  Prom 
the  numerical  results  represented  by  Fig.  2,  it  is  expected  that  both  XL  ^  and  XL  2 
approach  zero  as  \  approaches  zero.  Hence,  XL  ^  and  XL  2  have,  an  expansion  of 
the  form: 

_  3 

XL  —  C  ^  ^  ^  t  ....  (7l) 

where  only  odd  powers  occur  because  the  XL.  ^  are  odd  functions  of  \  .  The  values 
of  c  and  d  are  determined  by  substituting  this  expansion  in  eq.  (25)  and  equating 
the  coefficient  of  each  power  of  ^  appearing  to  zero. 


The  terms  ~  and  A  in  the  transcendental  equation  (25)  are  given  in  full  by 

eqs.  (23)  and  (2U).  Expanding  the  Bessel  functions,  & 

may  be  written  in  the  form; 
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In  order  to  be  able  to  clear  eq.  (25)  of  fractions,  the  quotient  of  the  two 
infinite  series  in  eq.  (72)  is  written  as  a  single  series: 


A  = 


l 


Using  this  value  the  transcendental  eq.  (25)  may  be  put  in  the  form: 

•  •  •  •  J  -4- 
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Now  the  assumed  expansion  given  by  eq.  (71)  is  substituted  above  for  XL  •  The 
result  is  a  series  of  the  form: 

b,  f*+M6  +  'M8  +  ••••  =  °  (75) 

Each  coefficient  b  of  this  series  must  be  zero  since  the  value  XL  is  a  root 

n  h2 

of  the  transcendental  eq.  (25).  Equating  b^  to  zero  and  neglecting  compared 

to  unity  yields  an  equation  involving  c  but  not  d.  This  is  interpreted  as  an 
equation  on  the  yet  unknown  c.  The  roots  are: 


Ct  i  c  t  C- 


where 


2AR  -V  R  +  Rx(\-rz)  +  J  [l  AR  +  ^  +  R1  (l-^1)]*  -  4  R1  0  -V  r)(?  At  K) 
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Equating  bg  to  zero  yields  an  equation  involving  both  c  and  d.  This  equation 
is  considered  to  determine  two  values  of  d,  one  for  c  ■  c^  and  one  for  c  »  c^. 

The  solution  of  the  equation  b^  ■  0  for  d  yields: 


d.  = 


C  CL 


^(~)  (2  A  4T?)  +  8  R(2  A  +  l)  +-  8Rl('  -V7-) 


(78) 
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